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Abstract 

Background: The myocardium exhibits heterogeneous nature due to scarring after 
Myocardial Infarction (Ml). In Cardiac Magnetic Resonance (CMR) imaging, Late 
Gadolinium (LG) contrast agent enhances the intensity of scarred area in the 
myocardium. 

Methods: In this paper, we propose a probability mapping technique using Texture 
and Intensity features to describe heterogeneous nature of the scarred myocardium in 
Cardiac Magnetic Resonance (CMR) images after Myocardial Infarction (Ml). Scarred 
tissue and non-scarred tissue are represented with high and low probabilities, 
respectively. Intermediate values possibly indicate areas where the scarred and healthy 
tissues are interwoven. The probability map of scarred myocardium is calculated by 
using a probability function based on Bayes rule. Any set of features can be used in the 
probability function. 

Results: In the present study, we demonstrate the use of two different types of 
features. One is based on the mean intensity of pixel and the other on underlying 
texture information of the scarred and non-scarred myocardium. Examples of 
probability maps computed using the mean intensity of pixel and the underlying 
texture information are presented. We hypothesize that the probability mapping of 
myocardium offers alternate visualization, possibly showing the details with 
physiological significance difficult to detect visually in the original CMR image. 

Conclusion: The probability mapping obtained from the two features provides a way 
to define different cardiac segments which offer a way to identify areas in the 
myocardium of diagnostic importance (like core and border areas in scarred 
myocardium). 



Introduction 

Patients who have suffered but survived myocardial infarction (MI) may subsequently 
suffer a possibly disabling or fatal cardiac arrhythmia. Late Gadolinium (LG) enhanced 
Cardiac Magnetic Resonance imaging (CMR) is used for assessing morphology of the 
myocardium in patients after ML An example image is depicted in Figure 1 with 
myocardium and scar delineated by green and red lines, respectively. The myocardial 
muscle fibers are completely dead at the core area of the scar. Thus, the core area does 
not react to the electrical signals propagating through the heart muscle telling the heart 
to contract. 

O© 201 3 Kotu et al.; licensee BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative 
BiolVted Central Commons Attribution License (http://creativecommons.Org/licenses/by/2.0), which permits unrestricted use, distribution, and 
reproduction in any medium, provided the original work is properly cited. 
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Figure 1 Cropped short-axis CMR image showing manual segmentation of the myocardium and scar 
tissues. The yellow and blue dots in the image are manually marked (by Cardiologist) coordinates to 
segment myocardium and scar. The red and green contours generated by cubic spline interpolations of the 
above coordinates show myocardium and scar tissues respectively. 



At the border areas, however, the muscle fibers are not all dead. The electrical signals 
will be disturbed in these areas possibly causing reentries and sometimes arrhythmias. 
The latter gives us reason to believe that the border areas are a very important part 
of the scar, and thus good ways of defining and visualizing border areas would be 
beneficial in many situations. Previous works [1,2] emphasize the importance of the 
border area (also referred to as peri-infarct area, or gray zone area) being a major 
determinant for identifying patients at risk of arrhythmia. The heterogeneous nature of 
cardiac segments is related to various diagnoses and prognosis. For example, Risk of 
serious arrhythmias and better indication for implantation of implantable cardioverter 
defibrillator (ICD) [3], heart rate of ventricular tachycardia [4], risk and degree of heart 
failure, indication of implantation of pacemakers in order to synchronize cardiac con- 
tractions [5] etc. Our main focus in this work is to develop methods to visualize, 
define and analyze the myocardium to identify regions of diagnostic importance in the 
myocardium. 

Earlier reported work 

To our knowledge, visualization of the heterogeneous nature of the myocardium is mostly 
based on thresholding methods described in the literature where thresholds are defined at 
intensity levels corresponding to some percentage of the maximum intensity level inside 
the scar area [1-3]. Recently, Shokrollahi etal. was able to divide the scarred myocardium 
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of animals into core and gray zones using fuzzy clustering algorithm and validated the 
results using in-vivo electro-anatomical voltage mapping [6]. Some previous publications 
reported automatic segmentation of scarred myocardium in late enhanced magnetic res- 
onance images [7-9]. The segmentation of scar itself is an important application as the 
scar size contains information useful for finding the patients with high and low risk of 
arrhythmia [3]. Also, scar size is largely responsible in ventricular remodeling [10]. Dicki 
et al use support vector machine (SVM) to classify the scar and the myocardium in 
delayed enhancement magnetic resonance images in [7] using three features. The three 
features used in [7] depend on the intensity of scar, myocardium and the whole CMR 
image. 

From our groups previous work [11], it is shown that patients with high and low 
risk of getting arrhythmia can be classified using textural measures from gray-level 
co-occurrence matrices. In our recent work [12], we presented some preliminary results 
where dictionary learning and sparse representations were used for representing the 
texture in the myocardial segments followed by a classification. 

The proposed probability mapping 

All the above discussed previous work have sought to identify methods to segment 
the hyper-enhanced area of the myocardium. However, by only delineating the hyper- 
enhanced myocardium, important information concerning the tissue within the infarct 
territory is not provided. Hence, our focus is on developing a method that provides prob- 
ability values as a means to distinguish the areas of the myocardium rather than the crisp 
definitions offered by a scar vs no scar segmentation. 

In this work, we propose a technique for transforming the scarred myocardium in the 
CMR image to a probability map to indicate varying degrees of its functioning capabili- 
ties: where high and low probabilities will indicate the scarred tissue and normal tissue, 
respectively. Section 'Calculation of probability maps using Bayes rule' describes the gen- 
eral principles for calculation of the probability maps, where each pixel is the posterior 
probability of that pixel being scar compared to healthy myocardium. The posterior prob- 
abilities are calculated using Bayes rule with the help of features from the myocardium. 
Features are calculated pixel by pixel based on a local neighborhood and is explained in 
Section 'Extraction of features'. We have based our experiments on two features; a local 
mean (direct current, DC) that will correspond well to how the doctors manually segment 
the scar, and a textural feature using sparse representations, that might reveal information 
regarding gray zone areas. 

Section 'Experiments and results' shows the experimental framework and examples of 
probability maps for both features. Also in section 'Experiments and results', cardiac seg- 
ments example is demonstrated to show the diagnostic importance of texture feature 
probability maps. Finally, the paper is concluded in section 'Conclusion. 

Calculation of probability maps using Bayes rule 

The class specific probability density functions (PDFs) modeling specific feature vector 
values known to be located within scar or healthy myocardial areas are estimated as 
p(v\scar) and p(v\myo), respectively. According to Bayes rule [13], the probability of a 
pixel being in the scar area is related to the density function and the prior probability 
function: 
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(P(scar)p(v\scar)) 

P(scar\v) = — , (1) 

piy) 

where piy) = P(scar)p(v\scar) + P(myo)p(v\myo) ) and v is a feature vector. The prior of 
scar, P(scar), reflects the mean value of the scar size in the training data set, and Piscar) 
is calculated as the ratio of the number of pixels present in the scar with respect to the 
myocardium in the entire training data. 

The parameters in the class specific PDFs can be estimated using parametric or non- 
parametric methods, and here they are calculated using maximum likelihood (ML) 
estimation [14] which is a method of low complexity. 

Our primary focus is deriving the probability mapping, and there is a strong link to 
segmentation as we define the functions so that they should be optimal in the sense of 
discriminatively expressing the properties of the scar and normal tissue. The minimum 
error classifier corresponds to classifying as scar when P(scar\v) > P(myo\v). Therefore, 
we evaluate the features and feature combinations used in various mapping with respect 
to their segmentation capability. As we will discuss later, this might not be entirely fair as 
in our case; features reflecting the cardiologists crisp decisions are favored. 

The extraction of feature vectors for myocardium is illustrated in section 'Extraction of 
features'. Feature vector v of the myocardium in general can be represented as: 

v = [feature v feature^ . . . feature^ , (2) 

where d is the dimension of the feature vector v. Finally, the probability maps are obtained 
by regenerating a new image such that each pixel in the myocardium is given by the value 
of P{scar\v) as calculated for that pixel. The probability values range from 0-1 and they 
are displayed using a color code shown in Figure 2(b) and (c). The shades of red and 
yellow in the color code represents high probability of scar and shades of green and blue 
represents low probability. 

Training and testing 

The process of generating probability maps is depicted in Figure 3. The whole CMR data 
set is divided into training and testing images. In the training process, training features 




Figure 2 Probability maps of the scarred myocardium in CMR images of a patient implanted with ICD 
using the DC and texture features, (a) The cropped CMR image containing left ventricle, (b) Probability 
map using the DC value as feature: Dark Red 1 00%-90%, red 90%-80%, yellow 70%-60%, green 60%-50%, 
cyan 50%-30% and blue 20%-0% probability of being scar. The other image parts are shown in light pink 
color, (c) Probability map using the texture feature. The same color code described above is used here. White 
contour shows the cardiologists segmentation of scar. 
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Figure 3 Training (dashed lines) and test (solid lines) phases involved in computing the probability 
maps of the scarred myocardium. 



are extracted along with some training parameters required in the testing phase. In the 
next step, the calculated training feature vectors are used for estimating the PDFs of scar 
and myocardium region using ML. The training process and its outputs are indicated by 
the dashed lines. The probability maps are generated on test images, different from the 
training images, and this process is indicated by solid lines. The extraction of features for 
training and testing phases are explained in section 'Extraction of features'. 

Extraction of features 

As discussed in section Introduction, two features, DC value and texture feature are 
used to generate probability maps of the scarred myocardium in this work. The detailed 
process of extracting the two features are discussed in the following sections 'DC Feature' 
and 'Dictionary-based textural feature'. 

DC Feature 

Historically, in electronics field, the mean is commonly referred to as DC (direct current) 
value [15]. Because of the intensity differences between the healthy and scarred 
myocardium tissues in the LG enhanced CMR images, DC value dc(i,j) is one of the fea- 
tures explored for developing the probability maps. This is probably the feature that would 
mimic the cardiologists segmentation best since they use the intensity values of scarred 
myocardium in CMR image to segment the scarred areas manually. We define the fea- 
ture dc(ij) = mean(I ^ x ^(hf)) (where / is the neighborhood around the pixel 
/(/,/)), such that a new image 7^ c , with dc(i,j) as values at pixel position is made with 
a sliding window averaging over the image. 

Dictionary-based textural feature 

A dictionary D, is an ensemble of a finite number of atoms, d^ k \ and can be used to 
represent signals. A linear combination of some of the atoms in the dictionary gives exact 



Kotu etal. BioMedical Engineering OnLine 201 3, 12:91 
http://www.biomedical-engineering-online.eom/content/1 2/1 /91 



Page 6 of 19 



or approximate representation of the original signal Let a column vector x of finite length 
N represent the original signal The dictionary atoms are arranged as columns in an N x K 
matrix D where K > N. The representation or the approximation of the signal, x and the 
representation error, r can be expressed as: 

K 

x = w 0<)d (k) = Dw, r=x-x=x- Dw, (3) 
k=i 

where w is sparse coefficient vector. In a sparse representation, only a small number of 
the coefficients w(k) are allowed to be nonzero. Finding the sparse coefficient vector can 
be a NP-hard (Non-deterministic Polynomial-time hard) problem and formulated as: 

w opt = argmin\\x - Dw\\ 2 s.t \\w\\q < s. (4) 

w 

The pseudo-norm || . ||o is the number of nonzero elements and, s denotes the sparseness. 
A vector selection algorithm called Order Recursive Matching Pursuit (ORMP) algorithm 
is used for sparse approximation [16]. 

Dictionary learning is the task of learning or training a dictionary on an available train- 
ing set such that it adapts well to represent that specific class of signals. The training 
vectors, L and the sparse coefficients are arranged as columns in the matrices X and W, 
respectively. The objective in dictionary learning is to give a sparse representation of the 
training set X in order to minimize the sum of the squared error. The cost function is 
formulated as: 

argmin F(D, W) = argmin \\ri\\^ s.t ||w|lo < s. (5) 

d,w d,w • 

This is a very hard optimization problem. Recursive Least Squares Dictionary Learning 
Algorithm (RLS-DLA), [17] is used for dictionary learning in all our experiments. A brief 
explanation of the algorithm is given in Appendix. 

Sparse representations and learned dictionaries have been shown to work well for tex- 
ture classification by Skretting and Husoy in [18] and by Mairal et al. [19]. Texture in a 
small image patch can be modeled as sparse linear combination of dictionary atoms in 
Frame Texture Classification Method (FTCM) presented by Skretting etal. in [18]. FTCM 
is developed by modeling a texture as a tiled floor where all the tiles are identical The 
color or the gray level, at a given position in the floor is given by an underlying continuous 
periodic two-dimensional function. It is shown based on this model that a vector from 
spatial neighborhood is indeed a sparse combination of finite dictionary atoms. 

The algorithm for sparse representation in our work proceeds as follows in the train- 
ing phase: Consider the myocardium in a CMR image, / that contains two texture classes: 
healthy and scarred myocardium. The training vector, yi of length N is made from that 
specific pixel and its neighborhood, -s/N x -s/N for each pixel in the training image. In the 
training set, each pixel is classified into a specified texture class. Then, the dictionaries, 
D s and D m are trained for the predefined texture classes (scar and myocardium) using 
RLS-DLA. Every training vector, yi is then represented sparsely using ORMP vector selec- 
tion algorithm [16] twice, both with D s and D m . The residuals (or representation errors), 
R s and R m , are calculated for each pixel in the myocardium of a training image as: 

Rsihj) = Wyi ~ A"? II and Rmihj) = Wyi - D m w™\\, (6) 
where zxv&wj 1 are sparse coefficient vectors. 
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The residuals, R s and R m , are combined to form the texture feature used in this work. 
The texture feature is calculated as the ratio of the residual of D m to the sum of the 
residuals of D s and D m : 

Rpihj) = Rm{Uj)l{Rs{i,j) + R m (i,j)). (?) 

Pixel by pixel R p value can be interpreted as the scaling of residuals R s and R m . Smaller 
values of R p means that the pixel is not likely to be scar (i.e. healthy myocardium), and 
larger values means that it is likely to be scar or border area. The texture feature, R p (i,j) 
from the training set is used to estimate the PDFs for the probability maps computation. 
In the testing phase, texture feature vectors are collected from the test image set, in the 
same way as training feature vectors. In Figure 3, training and classification parameters 
are indicated. In fact, these training parameters are the dictionaries, D s and D m . Using D s 
and D m learned in the training phase, the residual images, R s and R m are calculated for 
test images. Texture is not pixel-by-pixel local except in the edge between two textures; 
thus some sort of smoothing is used on the test residual images before calculating the 
probability maps [18]. The test residual images R s , and R m are smoothed using an a x a 
pixels separable Gaussian low pass filter with variance a 2 before calculating R p . 

Experiments and results 

The CMR images used in our experiment were provided by the Department of Cardiology 
in Stavanger University Hospital and were obtained from 44 patients with myocardial 
infarction. 24 patients had old MI and ICD was implanted in all patients as primary- or 
secondary prophylaxis. This group was defined as the high risk arrhythmic group. The 
remaining 20 patients had MI one year prior to CMR imaging and during this observation 
period, no reported incidents of serious arrhythmias occurred. This group was defined as 
the low risk arrhythmic group. The study was approved by the Regional Ethics Committee 
and informed consent was obtained from all patients. 

All the CMR images were obtained from 1.5 Tesla Philips Intera machine. Images were 
acquired with a pixel size of 0.82 x 0.82 mm 2 , covering the whole ventricle with short- 
axis slices of 10 m thickness, without inter-slice gap. The LG enhanced CMR images were 
stored according to the digital imaging and communications in medicine (DICOM) for- 
mat with 512 x 512 pixel resolution. The number of image slices with visible scar in each 
patient varies depending on the size of the scar and the scar size varies from one slice to 
the other. Only short- axis image slices were used in our experiments. Out of the total 444 
CMR images belonging to the 44 patients, 11 CMR images of poor quality due to motion 
and artifacts were not included. An experienced cardiologist delineated the borders of 
the myocardium and the infarction tissues as shown in Figure 1. The resulting set of seg- 
mented images were further split into labeled training and test sets. All our experiments 
were carried out in MATLAB. The CMR training and test images were cropped from 
512 x 512 original CMR image so that they contain only left ventricle with myocardium 
and infarction areas to save execution time. Preprocessing of any kind was not used on 
these cropped images. 

In all the CMR images, we took into account only the myocardium segmented by the 
cardiologists. Patients with high risk of getting arrhythmia have profound scars with pos- 
sibly large gray-zone areas. We therefore anticipated that images from these patients 
contained the information we wanted to represent in the dictionaries: information from 
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the myocardium, the core area of scars, and the gray-zone areas. Hence, the high risk 
patients group was used in the training phase. From the group of 24 high risk patients, 
6 patients were randomly selected to form a training set, whereas the remaining 18 high 
risk patients were used to form a test set. Two sets of training vectors were generated, one 
from the scar, and one from the healthy myocardium. 

This section is organized as follows: Sections 'Probability maps obtained using DC 
feature, dc(i,jY and 'Probability maps obtained using texture feature, R p (i,j)' explains 
the processes of calculating the probability maps using DC and texture feature respec- 
tively. Section 'Examples of probability maps obtained from DC and texture feature' 
discusses the examples of probability maps obtained in our experiments. Section 
'A possible application - cardiac segments' illustrates how probability maps obtained with 
texture feature can be used to define cardiac segments. Section 'The discriminatory power 
of the DC and texture features' considers the discriminative power of the DC and texture 
features. 

Probability maps obtained using DC feature, dc(i,j) 

The neighborhood size 3x3 was used to form training vectors as explained in 
'DC Feature'. The same neighborhood size must be used while training and finding the 
DC images 7^ c . The DC images obtained from the training images were used to form the 
training feature set. The parameters of the class specific PDFs: the mean and the stan- 
dard deviation were found using ML estimation using the training feature set. DC- values 
were scaled to have zero mean and unit variance before finding the ML estimates. The 
scaling coefficients from the training were stored to scale the test vectors. The probability 
maps of the scarred myocardium were calculated using Bayes rule as explained in section 
'Calculation of probability maps using Bayes rule'. An example of a probability map 
obtained using the DC feature with color code is shown in Figure 2(b). Note the sharp 
transition from the segmented scar to the myocardium. 

Probability maps obtained using texture feature, R p (i,j) 

The training vector and test vectors were generated in the same way as in the DC feature 
experiment using the same neighborhood size 3x3. Consider the pixels on the border 
zones, their neighborhood extends into other regions that are not under consideration. 
If we use training vectors from border regions, the dictionaries might learn the texture 
properties of other regions along with the texture properties they are intended to learn. 
So, the training vectors for the pixels whose neighborhood span other regions were not 
considered in our experiments. This is depicted in Figure 4. 

The dictionaries were learned on the training vectors generated from scar and healthy 
myocardium using RLS-DLA. The dictionary size of 9 x 90 atoms was used in our exper- 
iments. The initial dictionaries were formed by randomly selecting 90 vectors of length 9 
from the training sets. Dictionary size was selected based on previous experience, and it 
should be K > N << L. The forgetting factor is initialized to 0.995 and slowly increased 
towards 1 according to the exponential method described in [17]. If the number of dic- 
tionary atoms used to represent the image patch increases, then the residual decreases, 
and a large number of atoms will provide a good approximation with any full rank dic- 
tionary. Therefore, higher sparsity lowers the difference between the residuals of healthy 
myocardium and scar which in turn decreases the differentiation between the classes. 
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Figure 4 Training vectors are generated from each pixel as long as that entire neighborhood lies 
within one texture region. The neighborhoods of pixels P] , Pj and P3 include more than one texture, and 
the corresponding feature vectors are excluded from the training set. P4 and P5 have the entire neighborhood 
within one texture, and hence the corresponding feature vector is included in that texture's training set. 



The sparsity s used in our work is two, i.e. we used two vectors from the dictionary to 
represent the original image patch. 

After learning the dictionaries, they were used to obtain the residual images. For each 
pixel in the myocardium of training images, the scaled residual R p (i,j) was found from 
residual images as explained in section 'Dictionary-based textural feature'. The training 
feature set was then generated from these scaled residual values. The parameters of the 
PDFs to be used with the classifier was estimated from the labeled training vector set using 
ML estimation. The test vector set in this experiment was calculated after smoothing 
(using low pass Gaussian filter with a = 5 and 9x9 window size) the test residual images. 
Using the ML estimates and Bayes rule, the probability maps for texture features of the 
scarred myocardium were obtained. An example of such a probability map with color 
code is shown in Figure 2(c). The color code used is explained under Figure 2. Here, the 
dictionaries were learned without the removal of the DC value in each image patch. 

Dictionary learning was performed with and without the removal of the DC factor from 
the image patch of the training set. In order to illustrate some of the texture captured 
by the dictionary atoms, they were plotted as image patches in Figure 5. The difference 
between the scar and the healthy myocardium dictionaries was visible when the DC value 
was not removed, and this is according to the known intensity difference between scar 
and healthy myocardium. To investigate the variation in the dictionary atoms used as 
image patches, we had calculated a measure from the first derivative of dictionary atoms 
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Figure 5 Dictionaries of the healthy myocardium and scar: D m and D s with DC value and removing DC 
value of 3x3 patches. Absolute mean of the first derivative of the dictionary atoms: a) D m with DC- 0.3972 
b) D s with DC- 0.1455 c) D m without DC- 0.7314 d) D s without DC- 0.7336. 



in both vertical and horizontal directions. The absolute mean of the first derivative of 
the dictionary atoms was used to compare the scar and healthy myocardium dictionaries 
learned with and without DC value. The calculated absolute means are shown under 
Figure 5. The difference between the scar and healthy dictionaries without DC removal 
was observed clearly, both from the depicted dictionaries and the absolute means. The 
difference between the scar and healthy myocardium dictionaries trained by removing 
the DC value from image patches could not be observed visually. A larger block for 
DC removal prior to textural capturing is probably more appropriate, more similar to 
suppression of difference in illuminance conditions in computer vision. 

Figure 6 shows probability maps of some example CMR images and all of them were 
calculated using the texture feature, R p (i,j). In the middle column, the DC value was not 
removed from each image patch before training dictionaries and calculating the probabil- 
ities. In the last column, the DC value was removed from every image patch both before 
learning dictionaries and before finding probabilities. The DC value is obviously impor- 
tant, so we do not expect this experiment to provide competitive results, but we see a 
tendency of higher values in the scar and border areas, but not necessarily in the core of 
the scar. We expect the scar core to be very bright and could be captured easily by some 
sort of intensity based measure (for example the DC feature). If the physical nature of the 
tissue is reflected in the image texture, we can expect the scar core to be quite homoge- 
neous. Thus, we are not surprised that the regions that are the brightest in column II were 
actually quite dark in column III (core), but the intermediate regions were captured fairly. 
We have continued with the experiments in this work by letting the dictionaries incorpo- 
rate the DC value as we did not combine the DC feature, dc and texture feature, R p into a 
2D dimensional feature for computing the probability maps. Unless it is mentioned under 
figures of probability maps, all the probability maps of texture feature, R p were computed 
with dictionaries learned without removing DC from every image patch. 

Examples of probability maps obtained from DC and texture feature 

The second and the third column of Figure 7, shows more examples of probability maps 
of scarred myocardium using the DC, dc(i,j), and texture, R p (i,j), features, respectively. 
Figure 8 shows examples of probability maps of scarred myocardium for both DC and 
texture feature for the patients without ICD. This shows that the PDFs estimated for the 
DC and texture features from the patient group implanted with ICD were also able to 
produce probability maps for low risk arrhythmia patients (i.e. patients without ICD). 
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Figure 6 Probability maps of the scarred myocardium in CMR images using the texture feature, R p . 

Column I: The cropped CMR images of 3 patients implanted with ICD containing left ventricle. II: Probability 
maps of the texture feature computed with dictionaries trained on image patches (used for collecting 
training vectors) without removing the DC value. Ill: Probability maps of the texture feature computed with 
dictionaries trained on image patches by removing the DC value. White and black contours shows the 
cardiologists segmentation of scar. 



Both the features were able to identify the healthy myocardium without scars, and an 
example of this is shown in Figure 9. In order to compare the probability maps obtained 
by the two methods, they are plotted with the same color map. Experiments of probability 
mapping combining the DC feature and texture feature gave almost identical probability 
maps as the DC feature because of the dominance of the DC feature over the texture 
feature. Improved ways of combining both the features will be explored in our future 
work. 

A possible application - cardiac segments 

A cardiac segment, CS can be defined as a region of interest in the myocardium which 
might have diagnostic importance. We define a candidate cardiac segment by a lower and 
upper probability so that the segment consists of all pixels with probability values in the 
range between these. Different candidate cardiac segments might be evaluated, for exam- 
ple, calculating one or several features from the candidate segment and correlate result to 
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I II III 



Figure 7 Probability maps of the scarred myocardium in CMR images of 5 patients implanted with 
ICD using the DC and texture features. Column I: The cropped CMR images containing left ventricle. II: 
Probability maps using DC value as the feature. Ill: Probability maps using the texture feature. Refer to the 
Figure 2 for the color description. 
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Figure 8 Probability maps of the scarred myocardium in CMR images of patients without an ICD 
using the DC and texture features. Column I: The cropped CMR image containing left ventricle. II: 
Probability map using DC value as the feature. Ill: Probability map using the texture feature. Refer to the 
Figure 2 for the color description. 



a clinical meaningful hypothesis. We demonstrate the concept by making a comprehen- 
sive search for significant cardiac segments by comparing two patient groups and testing 
for statistical significance. The demonstration experiment is meant to illustrate how the 
probability mapping can be used to localize cardiac segments containing information 
discriminating between the high and low risk patient groups. 

For each patient, cardiac segments were determined by a lower boundary L and an 
upper boundary U. These boundaries were varied in the range 0-1 in steps of 0.025 with 
the restriction that U — L > 0.1. Figure 10(a) illustrates this by visualizing each of the 
resulting candidate cardiac segment as a black dot with coordinate (X, U). For example, 
the three yellow dots in coordinates (0.15,0.3), (0.3,0.7), and (0.7,0.825) correspond to 




Figure 9 Probability maps of the myocardium without scar in CMR images of patients with ICD using 
the DC and texture features. Column I: The cropped CMR image containing left ventricle. II: Probability 
map of the myocardium without scar using DC value as the feature. Ill: Probability map of the myocardium 
without scar using the texture feature. Refer to the Figure 2 for the color description. 
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Figure 10 Exploring cardiac segments of a patient implanted with ICD. (a) The coordinates (L, U) shows 
all cardiac segments (CS) corresponding to the intervals L — U. (b) An image slice from a patient from the 
group with high risk of arrhythmia, (c) The corresponding probability map, and (d-f) the cardiac segments 
corresponding to the yellow (/_, U) coordinates CS] , CS 2 , and CS3. Subplots (g-i) shows the corresponding 
cardiac segments derived without probability mapping. 



the cardiac segments CS\ = 0.15 — 0.3, CS2 = 0.3 — 0.7, and GS3 = 0.7 — 0.825, respec- 
tively. Figure 10(b) and (c) show one of the original image slices and the corresponding 
probability mapping (texture feature) for one of the high risk patients. The probability 
maps, which do not have 0 — 1 probability range, was extended using a sigmoid function 
for easy comparability. Figure 10(d-f) shows the cardiac segments corresponding to CS\, 
CS2 and GS3 for the image slice shown in Figure 10(b). 

Preliminary experiments indicate that corresponding segments determined without 
using the probability mappings depends on the signal intensity values and seems to lack 
the ability to define contiguous areas in the same way as the probability mapping for 
some of the cardiac segments. For example, the cardiac segment determined directly from 
thresholding relative to the maximum 3D scar signal intensity is shown in Figure 10(g-i) 
for boundaries corresponding to CSi CS2 and GS3, respectively. These plots are com- 
parable to the corresponding cardiac segments based on the probability mapping in 
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Figure 10(d-f), and show less contiguous cardiac segment, especially CS\. Figure 11 
shows similar plots of the cardiac segments corresponding to CS\, CS2 and GS3 for one 
of the low risk patients. 

For each candidate cardiac segment, its relative size (accumulated for all slices) to 
the total myocardium was calculated. The computed relative sizes for high and low risk 
patient groups were compared using a Mann-Whithney test at significance level 0.05. The 
LU-plots in Figure 12 show the p-values for the test of all the candidate cardiac segments 
derived using the texture based (a) and DC based (b) probability map and those derived 
without the use of probability mapping (c). The blue colored dots indicate cardiac seg- 
ments for which the high and low risk groups showed significant differences of the relative 
size. 

Interpreting the comparison results for the texture feature probability map based car- 
diac segments in Figure 12 (a), the significant CS with L = 0 and U > 0.5 correspond 




Figure 1 1 Similar plots and cardiac segments as in Figure 1 0 for a patient from the group of patients 
with low risk of getting arrhythmia, (a) The coordinates (L, U) shows all cardiac segments (CS) 
corresponding to the intervals L — U. (b) An image slice from a patient from the group with low risk of 
arrhythmia, (c) The corresponding probability map, and (d-f) the cardiac segments corresponding to the 
yellow (/_, U) coordinates CS] , CS 2 , and CS 3 . Subplots (g-i) shows the corresponding cardiac segments 
derived without probability mapping. 
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Figure 1 2 Plot of (L, U) coordinates showing all cardiac segments (CS) corresponding to the intervals 
L - U comparing the relative sizes of the cardiac segments for two patient groups. The color of the 
(L, U) coordinates indicates the p-value resulting from using the Mann-Whithney statistical test to determine 
if the relative sizes of the cardiac segments are different in the two patient groups. Subplots (a) and (b) 
corresponds to the tests of cardiac segments derived with probability mapping of texture and Defeatures, 
respectively. Subplot (c) shows the tests of cardiac segments derived directly from the original image. 



to the cardiac segment containing all of the healthy myocardium and gradually including 
the scar area. More interesting is the localization in the upper right corner which makes 
more specific localization in the myocardial tissue in the regions close to the scar. Note 
that GS3 is identified as significant while CS\ and CS2 are not. Figure 12 (b) shows the LU- 
plot from the probability maps of DC feature and, the plot shows that cardiac segments 
does not give any localization. This is due to the narrow transition between high and low 
probabilities of the DC based probability maps. The DC probability maps will predomi- 
nantly have extreme values close to 0 or 1 and few values in between 0-1. This provides 
little information to explore. These plots (Figure 12 (a) and Figure 12 (b)) of cardiac seg- 
ments clearly shows that the probability maps obtained from the DC features reflect the 
way cardiologist perceive the scarred myocardium in CMR images and, the probability 
maps from texture features gives information about the characteristics of myocardial tis- 
sue and might emphasize information that is not easy for inspection with a human eye 
in the original CMR image. The difference between the narrow and broad transitions for 
the DC- and texture-based probability mappings might explain this. As the narrow tran- 
sition gives little information besides the crisp core/scar discrimination corresponding 
to the cardiologists delineation based on the visual information. Thus, most of the can- 
didate cardiac segments will be empty as there are virtually no intermediate probability 
values. This is not the case for the texture based probability mappings, which discloses 
the information present in this intermediate probability range. 

In comparison the results for the cardiac segments derived directly without using the 
probability plots (Figure 12 (c)), the localization seems vaguer, not specifically identifying 
any coherent region in the myocardium. It is our belief that probability mapping might be 
used to identify cardiac segments that might be used in future studies. 

The discriminatory power of the DC and texture features 

The probability maps calculated for the two features, DC and textural feature, give dif- 
ferent information. The discriminatory power of the two features were explored by 
comparing to the manual segmentation of scar and healthy myocardium to obtain 
Receiver Operating Characteristic (ROC) curves. The ROC curves calculated for the DC 
and texture features on the CMR images of the 18 test patients with high risk of getting 
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arrhythmia are shown in Figure 13. The area under the average ROC curves is 0.9052 and 
0.8428 for the DC and texture features, respectively. From ROC curves, it is observed that 
the discriminative power of the DC feature is comparatively higher than that of the tex- 
ture features. The sensitivity, and specificity reported in [7] are plotted as a single point in 
Figure 13. This reported result is not calculated on our CMR data. This single point lies 
below within 95% confidence interval of DC feature. This reported segmentation result 
performed well compared to our average DC and texture ROC curves as Dikici et al. [7] 
does not include all the CMRI slices of a patient for training and testing where as in our 
work we include all the CMRI slices i.e. the scar volume. 

Conclusion 

In this paper, a new technique has been proposed for enhanced visualization of scarred 
myocardium in LG enhanced CMR images using probability mapping. Using Bayes rule 
and ML estimation, we obtained the probability maps of scarred myocardium using DC 
and texture features from each pixel in the myocardium. The probability values obtained 
from the DC features can be used for determining the scar size, and they reflect the 
way cardiologist perceive the scarred myocardium in CMR images. The probability maps 
from texture features gives information about the characteristics of myocardial tissue and 
emphasize information that is not easy for inspection with a human eye in the original 
CMR image. Using textural features in the probability mapping provide us with a tool 
for defining cardiac segments as well as give us features that can be used to distinguish 




0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 



1 -Specificity 

Figure 1 3 The average ROC Curves with 95% confidence intervals for the DC and texture features on 
all the CMR slices of thel 8 test patients with high risk arrhythmia. The ROC curves are compared to the 
results reported by Dikici (etal.) [7]. Area under the average ROC curves for the DC and texture feature are 
0.9052 and 0.8428, respectively. 
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between different groups of patients (patients high and low risk of getting arrhythmia). 
Our belief that the probability maps calculated by the texture feature captures the het- 
erogeneous nature of scarred myocardium is justified with the example of the cardiac 
segments experiment. The two features tested in this paper gave valuable, but different, 
information. In future work, benefit of combining these two features as well as other 
features will be explored. 

Intellectual Property Rights Patent pending. The Patent Application Number is 
1121307.1. All Intellectual Property Rights are reserved by the University of Stavanger 
and Stavanger University Hospital. Prekubator TTO is responsible for IP management 
and commercialization. 



Appendix 

Recursive least squares dictionary learning algorithm 

The dictionary learning problem can be solved in a two steps algorithm. Step 1) Find 
sparse coefficient matrix, W, using vector selection algorithms, keeping dictionary, D, 
fixed. Step 2) Keeping W fixed, find D. The dictionary update step depends on the 
dictionary learning method we choose to use. D and C are initialized randomly. 

In this paper, an on-line dictionary learning algorithm, Recursive Least Squares 
Dictionary Learning Algorithm (RLS-DLA) presented in [17] is used for dictionary 
learning and ORMP vector selection algorithm presented in [16] is used for sparse 
representation. The RLS-DLA updates the dictionary with the arrival of each new train- 
ing vector. In deriving updating rules for RLS-DLA, X t = [x\,X2, of size S x t, 
Wt = [ w\, W2t • • • t Wt] of size K x t and Q = (W t W t ) are defined for the 'time step' t 
At each time step the dictionary D t -\ and the C matrix are updated so that they obey the 
least squares solution D t = (X t Wf ) ( W t Wj) . The matrix inversion lemma (Woodbury 
matrix identity) is applied on Q to get the following update rules: 

Q = Q_i — auu T , (8) 

D t = D t -! + ar t u T , (9) 

where u = C t -\w t and a = 1/(1 + wju), r t = Xt — D t -\w t is the approximation error. 

With the inclusion of an adaptive forgetting factor X t in RLS-DLA, the update 
equation 8 is changed to: 

Q = (A^Q-i) -auu T , (10) 

where u = (A^~ 1 Q_i)m/^. The remaining equations remain the same. Due to introduction 
of the forgetting factor, RLS-DLA has good converging properties and is less dependent 
on the initial dictionary. 
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